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Abstract 

We study the motion of a massive particle in a quenched random environment at zero 
temperature. The distribution of particle positions is investigated numerically and special focus 
is placed on the mean stopping distance and its fluctuations. We apply a scaling analysis in 
order to obtain analytical information about the distribution function. The model serves as a 
simple example of transport in a random medium. 
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distributions. 



1 Introduction 

Dynamics in random media constitute a set of phenomena widely studied in modern statistical physics 
and soft condensed matter. The problems in questions extend from the growth of interfaces to the 
diffusion of scalars in disordered materials In order to understand the motion of manifolds 

in random environments with focus on the pinning mechanism, thermal fluctuations are to a first 
approximation often disregarded, see e.g., [0-0. A particular simple case is the motion of a massive 
particle, i.e., a zero dimensional manifold, in a random medium subject to pinning. 

This model was recently analyzed by Stepanow and Schulz |]rU| . They studied the behavior of a 
Newtonian particle of mass m in a random environment in d-dimensions described by the equation 
of motion mx = F(x); here F is a quenched random force with a Gaussian distribution of strength 
a. In dimensions d larger than one, the large time results of Stepanow and Schulz for the first two 
cumulants are given by 

/ \ VnTn ^ . . /\\9\ ^0 ^ \ 

(x)^oc^^^vo and ((x " (xM )^ oc t . (1) 

Here vq denotes a unit vector along the direction of the initial velocity vq. In one dimension the 
results in Eq. are undefined owing to the diverging denominators. This behavior refiects the fact 
that in one dimension the force F{x) can always be derived from a potential. The energy is thus 



conserved in this model and the particle will simply oscillate back and forth between two potential 
barriers. 

In the present paper we extend the results of Stepanow and Schulz||To[| by including a friction 



term in the equation of motion and thus allowing for random pinning in the one dimensional case. 
The paper is organized as follows. In Section 2 we define the model and write down the appropriate 
equations of motion both on the level of the stochastic equation of motion and the associated Kramers 
equation for the time dependent distribution. In Section 3 we outline the numerical analysis applied 
to the equation of motion. In Section 4 we carry out a simple scaling analysis, identifying the relevant 
parameters and the scaling functions describing the two first cumulants of the stationary distribution. 
Section 5 deals with the numerical results and we close the paper with a summary and conclusion in 
Section 6. 



2 The Model 

In one dimension the motion of a particle of mass m moving in a viscous medium and subject to 
a random force field F{x), depending on the actual position x of the particle, is governed by the 
Langevin equation 

dj^ X dx 

m— + ^— = F{x)+r]{t) . (2) 
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Here 7 is the coefficient of friction and rjit) a thermal Gaussian white noise with correlations 



T]{t)ri{t') = 2D5{t — t'). The noise strength D is given by the fluctuation-dissipation theorem, 
D = 'jkBT] T is the temperature of the medium. The distribution of F is also chosen to be Gaussian 
with the correlations 

= and {F{x)F{x)) = aS{x - x') ; (3) 

here a is the force correlation strength. 

In the present context we consider the case where the random trapping dominates the physics and 
we shall thus disregard the thermal fluctuations relative to the random environment characterized 
by F, that is we set D = in Eq. this limit corresponds to zero temperature or to a very massive 
particle m — >■ 00. Breaking up Eq. (0) it can then also be written as the set of equations 

dv , 

_ = -^v + F{x) , (4) 

where we for convenience have set m = 1. Equivalently, we can also choose to discuss the problem 
in terms of an equation for the distribution function P{x,v,t) itself, namely the Kramers equation 



in the limit of zero diffusion coefficient iO 

dP{x,v,t) 



P{x,v,t). (6) 



dt 

As earlier stated, in the limit 7 = the model described by Eqs. (|5|)-(||) above reduces to the 



one discussed by Stepanow and Schulz in [TU|. However, here we allow for dissipation by having a 
non-zero 7, and entirely different physics will arise. 



3 Discretization and Numerical Accuracy 

In order to render the problem amenable to numerical analysis we modify the force correlations 
Eq. defined on the continuum by applying a finite width to the delta function (or equivalently, 
formulating the problem on a lattice) in such a manner that 

{F{x)F{x')) = (a/ao) 5n{x),n{^') ■ (7) 

The subscript n{x) of the Kronecker delta in Eq. ([^ is thus a counter of which cell in the lattice the 
position X refers to; more precisely, n{x) = {n E Z \ x E [nao] {n + l)ao]}- In the limit Oq — > we 
recover the full continuum description in Eq. (|]). 

The equations of motion in Eq. are solved numerically by means of the fourth order Runge- 
Kutta method (RK4) [|12|. First, by calling a random generator an actual realization of the force 
field is generated corresponding to a specific environment. In the next step we solve the equations 
of motion with the specific initial conditions that the particle initiates its path at x = at t = 
with initial velocity vq. A typical trajectory originating from a simulation is shown in Fig. ^ below. 
Finally, we perform an average over realizations of the force field, typically of the order of 10000 
realizations. 
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There are three different sources of uncertainty in the simulations: The finite time step r, the 
finite lattice distance ao, and the finite number of samples of the random force distribution. Moreover, 
they are not independent: Decreasing ao will require a smaller r in order to keep the error small, 
since a smaller lattice constant leads to a more erratic potential in which the numerical solution is 
more sensitive to the magnitude of r. Choosing a value for we fix r such that the stopping distance 
of the particle in a specific environment does not change more than a prescribed fraction (numerical 
precision) when decreasing r by a factor of say 10. Then we sample a large number of realizations of 
the environment such that the statistical errors are below the desired precision. The typical relative 
numerical error in the simulations is about 1%, and the typical relative statistical error is no more 
than 5%. 

4 Scaling Analysis 

Before we embark on a more detailed discussion of the numerical results it is illuminating to perform 
a simple scaling analysis. Rescaling space and time in an affine manner according to x ^ xb = x and 
t ^ at = i, where a and b are scale parameters to be determined, the velocity scales like v vb/ a = v 
and we obtain by insertion in Eq. the scaled equation of motion dv/dt = bF{x)F{x) / a? — (7/a)?). 
The force on the rescaled lattice is thus given by F{x) = bF{x)a'^, is also Gaussian and correlated 
according to (F(x)F(x')) = b^ /a^a5{x — x'). Concluding, we infer the following scaling of the 
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parameters: a ^ a = ab^/a^, 7 — 7 = 7/a, fo — ?7o = {b/a)vo, and Oq ^ Oq = baQ. Choosing the 
scale factors a and b in such a manner as to ehminate the dependence on the force strength a and 
the friction coefficient 7, i.e., a = 7 and b = 7^/^/0"^/^, we finally obtain the scaling form 

<t) = ^G, (7t, ^0(7/^)'/', ao(7V^)'/') , (8) 

where Gi is a random function depending on the actual force realization. The above scaling analysis 
can, of course, also be carried out on the basis of the Kramers equation Eq. (^. 

In the following we limit our discussion to the statistical properties of the final position of the 
particle. In this limit the dependence on time disappears and we obtain 

X ^ x{oo) = ^G (vo ii/af' , ao{iV<yf") , (9) 

where G is a new scaling function. The physical meaning of the first argument in the scaling 
function G is related to energy considerations. In one dimension the force F derives from a potential 
U according to U{x) = — F{x') dx' and, correspondingly, the mean square fluctuation of U is 
given by (f/^(x)) = dx' dx" {F{x')F{x")) = ax. Consequently, the typical potential energy of 
the disorder behaves as y/ax] this also follows from the fact that the potential performs a 'random 
walk' in x-space. Evaluating this energy at the trapped position of a free particle, i.e., Vo/'y, we 
obtain a typical upper limit for the pinning energy. Comparing this energy to the initial energy of a 
free particle, ~ Vq, we arrive at the following parameter characterizing the degree of disorder: 



i.e., a power of the first argument in the scahng function. Here small values of the dimensionless 
parameter /i corresponds to strong disorder, whereas large values of fi characterizes weak disorder. 



5 Analysis of Numerical Data 

In this section we study the various moments of the variable x = x(oo) as a function of the initial 
velocity vq and the disorder strength a. Numerical accuracy and the associated computation time 
puts a limit to the range of values for which we examine the scaling function G. We have in particular 
analyzed the weak disorder regime determined by 10 < /i < 100. Stronger disorder requires better 
routines or much longer computation time in order to acquire the same accuracy. 

5.1 The trajectory 

It is instructive first to consider the outcome of a typical simulation. In Fig. 1 we have thus shown 
a simulated trajectory with a superposed trajectory of the free flight, corresponding to the solution 
of Eq. (§ with F{x) = 0, i.e., 

a;free(^) = —(1 " exp(-7t)) , (11) 

7 

with x(0) = and v{0) = vq. 

The initial motion follows that of a free particle, however, as the kinetic energy degrades due 
to friction the particle will soon experience a (random) force large enough to reverse the motion. 
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In order to illustrate this effect we depict Fig. ^.a the potential energy U{x) (solid curve) and the 
kinetic energy of free flight (the dashed curve) a function of position. The point of intersection of the 
two curves is to a first approximation the point of first reflection of the particle. The actual motion 
of the particle is, of course, sensitive to the disorder before the point of reflection, perturbing the 
dashed curve on the figure. After the first reflection the velocity is reversed and the process reiterates. 
Owing to the gradual loss of kinetic energy the particle travels shorter and shorter distances and 
finally comes to rest in a local minimum of the potential. Figure |^.b is a magnification of the region 
around the minimum that finally traps the particle. 

In Fig. ^a we show the distribution of final positions and note the general trend that the particle 
travels shorter than for the corresponding free motion. This is illustrated by the vertical line on 
the figure, indicating the final position for a free particle with the same initial velocity. Despite a 
superficial resemblance the distribution depicted in Fig. ^.a is not Gaussian. A transparent way to 
exhibit this deviation from Gaussian behavior is to plot the logarithm of the distribution; a Gaussian 
distribution then yields a parabolic shape. This plot is shown in Fig. ^.b, where we have also included 
a fit to a parabola. It follows from the fit that the actual distribution possesses long tails. 

5.2 Mean distance traveled 

A central quantity is the total mean distance traveled by a particle, (x) p, as a function of the initial 
velocity vo and the force correlation strength cr. The notation {■)f here denotes an average over the 
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Gaussian force distributions. In Fig. ^ we have depicted (x) i? as a function of a for fixed damping 
7 = 5.0 and initial velocity fo = 10.0. Data from two different seeds, i.e., two series of realizations 
of the random force field, are shown in order to indicate the statistical uncertainty. 

The data are reasonably independent of the small distance cutoff ag, and thus applying the scaling 
form Eq. (j^) we infer 

(^)F = ^(G'(t;o(7/^)^/')>^ . (12) 

A good fit to the data in Fig. ^ is provided by the stretched exponential {x)f ~ aexp (— 6cr'^). 
Averaging over a set of data with different cutoffs oq in the range from 0.01 to 0.001 we obtain for 
the parameters a = 2.000±0.001, b = 0.058±0.001, and c = 0.34±0.01. The error bars are extracted 
from an appropriate sample of data. The data are thus consistent with the scaling form 

{x}f exp , (13) 



7 V ^0 

which has the correct limiting form x ~ fo/7 in the absence of disorder or large initial velocity vo- 
Also, in the limit of vanishing friction 7—^0 the first moment {x)f vanishes. Consequently, there 
is a largest mean distance traveled for an appropriately selected value of the friction. Since we are 
here considering the weak disorder regime one could argue that the first two terms in the expansion 
of the stretched exponential are sufficient. We shall however stick to the expression in Eq. ([TB|) since 
it behaves correctly in the limit 7-^0. 



As an alternative check of the scaling form Eq. (13) we next vary the initial velocity Vq. The data 
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are shown in Fig. |^ for the values 7 = 5.0 and cr = 1.0. It follows that Eq. ( ]TB| ) does describe the 
data well. Also shown is a linear fit and a fit to a stretched exponential. We note that all curves 
appear to coincide on the figure, thus not allowing us to distinguish the best fit. 

5.3 Fluctuations in the stopping distance 

The last issue we address is the root mean square fluctuations of the stopping distance Ax = 
a/((x — depicted in Fig. |^. Here the situation is a little complicated since we find a dependence 

on the small distance cutoff qq. From the general scaling analysis yielding Eq. (P) we infer a scaling 
form for Ax, 

Ax = = [vo{llaf'\ a,{^V^Y") , (14) 

where G2 is a new unknown scaling function. 

In Fig. ^ we have for fixed damping 7 = 5.0 and initial velocity vq = 10.0 plotted the standard 
deviation of x versus the disorder strength a for two different cutoffs qq on a log-log scale. The linear 
fit suggests the simple power law behavior 

Ax = uia^\ (15) 

In order to determine the dependence on the lattice constant ao we have in Figs. |^ and ^ plotted 
the exponent zi and the amplitude ui as a function of ao in order to study the limiting function 
6*2(^0(7/0")^^^) = \imao^oG2{vo{'j/aY^^,aQ{'y^/ay^^) for decreasing oq. With a linear fit we find the 
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following values for the amplitude and the exponent: 

Ml = 0.143 ±0.001, ^1 = 0.326 ± 0.002, (16) 

indicating that the best fit is obtained using a constant scaling function G2 = const.. We are thus 
led to the form 

Ax = const.— , (17) 

independent of initial velocity. Substantiating our result we have finally in Fig. |^ plotted Ax versus 
vq. The independence of Ax on vq is quite convincing. 



6 Summary and Conclusions 

In this paper we have studied a simple model for damped ballistic transport in a random environ- 
ment: A massive particle moving in a random medium also subject to a viscous force. From energy 
considerations alone it is clear that the particle will eventually become pinned at a local minimum 
of the potential and we have thus restricted our attention to the properties of the stationary distri- 
bution. We performed an elementary scaling analysis of the equation of motion and exploited the 
results in the ensuing numerical analysis. Here we restricted our attention to the behavior of the first 
two cumulants, i.e. the mean stopping distance and the mean square fiuctuations of the stopping 
distance as functions of disorder strength and initial velocity. We found that the data in the weak 
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disorder regime were well described by the expressions: 

{x)f exp , {{x- {x)) ) ~ const.— . 18 

We, moreover, exhibited the non-Gaussian nature of the stationary distribution of stopping positions 
of the particle. 

There are many more quantities of interest that one could extract numerically. First of all it 
would be interesting to explore also the dynamic properties, that is the behavior of e.g. x(t) and v{t). 
Moreover, a more thorough analytical examination would also be desirable. Finally, an extension of 
the model to the case of a random force with a non-vanishing mean value would introduce some new 
facets into the problem. 
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Figure 1: We show a typical trajectory of ballistic motion in one dimension in a random 
force field (the full line). For comparison we also depict the trajectory of the motion in 
the absence of disorder (the dashed hne). Notice how the trajectories are very similar at 
short times, whereas for larger times they differ substantially 
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Figure 2: (a): Comparison of the typical energy of the disorder, which is a random 
walk in x-space, with the energy of the free particle, as a function of position, (b): The 
particle reflects a number of times on the potential and gradually looses its energy, until 
it is trapped in a local minimum of the potential. 
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Figure 3: In (a), we plot the stationary distribution of the particles rest position. The 
vertical line at x — 2 marks the final position of a free particle having the same initial 
velocity. In (b) we have shown the same data on a semi log scale. A fit to a parabola is 
also shown, and the deviation from Gaussian behavior is seen to be large especially away 
from the center. 
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Figure 4: The mean stopping distance as a function of disorder strength a. The data 
(from both seeds) have been fitted to a stretched exponential (sohd hne). 
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Figure 5: This is a plot of (x) versus Vq. We have shown both a hnear fit, a fit to 
stretched exponential, and also the prediction Eq. (0). They are all indistinguishable 
this plot. 
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Figure 6: A graph of standard deviation Ax of stopping distribution as a function of 
strength of disorder. The two different sets of data are for two different lattice constants. 
The straight hues indicate power-law behavior, albeit with different exponents and am- 
phtudes. 
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Figure 7: Exponent Zi (see Eq. (|T5| ) of tlie power law fits as a function of tlie lattice 
constant, and a fit to a straight line. We are interested in the limit Oq — > 0. 
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Figure 8: Amplitude Ui in Eq. (|^) of tlie power law fits as a function of the lattice 
constant, and a fit to a straight hne. 
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Figure 9: The dependence of the width of stopping distances on the initial speed Vq is seen 
to be well approximated by a constant expression. Prom this figure we have Ax — 0.154. 
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